1. Выбор данных

В качестве данных возьмём один из стандартных наборов для языка R, который описывает характеристики и оценки курсов и характеристики преподавателей Техасского университета.

2. Считывание и просмотр

Загрузим данные из файла и выведем первые 6 строк.

data <- read.csv("TeachingRatings.csv");
head(data)

3. Описание данных

В этом наборе данных есть 13 переменных:

  1. X — идентификатор курса;
  2. minority — принадлежит ли преподаватель к меньшинству;
  3. age — возраст преподавателя;
  4. gender — пол преподавателя;
  5. credits — является ли предмет факультативным;
  6. beauty — оценка физического облика преподавателя группой из шести студентов, усредненная по шести участникам оценки, нормированная;
  7. eval — оценка преподавания курса от 1 (очень плохо) до 5 (отлично);
  8. division — является ли этот курс курсом высшего или низшего отделения (Курсы низшего отделения — это в основном большие курсы первокурсников и второкурсников);
  9. native — является ли преподаватель носителем английского языка;
  10. tenure — имеет ли преподаватель учёную степень;
  11. students — количетсво студентов, участвовавших в оценке;
  12. allstudents — общее количетсво студентов, обучавшихся на этом курсе;
  13. prof — идентификатор преподавателя;

4. Типы признаков

  1. X — качественный;
  2. minority — качественный;
  3. age — количественный;
  4. gender — качественный;
  5. credits — качественный;
  6. beauty — порядковый;
  7. eval — порядковый;
  8. division — качественный;
  9. native — качественный;
  10. tenure — качественный;
  11. students — количественный;
  12. allstudents — количественный;
  13. prof — качественный;

Рассчитаем моду для признака “students”:

# Create the function.
getmode <- function(v) {
   uniqv <- unique(v)
   uniqv[which.max(tabulate(match(v, uniqv)))]
}
getmode(data$students)
## [1] 12

И посчитаем частоту моды:

length(data$students[data$students==12])
## [1] 21

Проведём аналогичные рассчёты для признака “allstudents”.

Мода:

getmode(data$allstudents)
## [1] 15

И частота:

length(data$students[data$students==12])
## [1] 21

Проведём аналогичные рассчёты для признака “age”.

Мода:

getmode(data$age)
## [1] 52

И частота:

length(data$students[data$age==52])
## [1] 47

Так как у признаков “allstudents” и “students” мода достигается 21 раз, а у признака “age” — 47, признаки являются дискретными.

5. Порядковые признаки

В этом датасете нет текстовых меток порядкового признака.

6. Matrix plot, outliers, etc.

pairs(~age+beauty+students+eval+allstudents,data=data)

## 7. Симметричность распределений. Сильно несимметричных (с хвостом вправо) распределений здесь нет.

8. Аутлаеры

По графикам видно, что, во-первых, allstudents и students имеют линейную зависимость (причём без outliers). Во вторых видно, что age и beauty имеют линейную зависимость, но менее явную.

9. Неоднородности.

По графикам видно, что где есть признак allstudents с другими признаками даёт сильную неоднородность.

Выведем те записи, которые им соответсвуют и попытаемся понять причину

data[data$allstudents>320, ]

В выборке есть очень популярный преподаватель

data[data$prof==73, ]

Уберём эти записи из выборки и посмотрим на скаттерплот

filteredData = data[data$prof!=73, ]
pairs(~age+beauty+students+eval+allstudents,data=filteredData)

#filteredData = data[data$prof!=73, ]
filteredData_s = data
filteredData_s['log'] = log(data$allstudents)
filteredData_s['log2'] = log(data$students)
pairs(~age+beauty+log2+eval+log,data=filteredData_s)

По графикам видно, что корреляция между количеством студентов, принявших участие в опросе и количеством студентов, прошедших обучение по курсу положительна и близка к единице.

Может большое количество студентов зависит от отделения?

pairs(~age+beauty+students+eval+allstudents,data=filteredData_s[filteredData$division=="lower", ])

pairs(~age+beauty+students+eval+allstudents,data=filteredData[filteredData$division=="upper", ])

Зависимость есть, но не особо сильная. А от учёной степени?

pairs(~age+beauty+students+eval+allstudents,data=filteredData[filteredData$tenure=="no", ])

pairs(~age+beauty+students+eval+allstudents,data=filteredData[filteredData$tenure=="yes", ])

pairs(~age+beauty+students+eval+allstudents,data=filteredData[filteredData[filteredData$division=="lower", ]$tenure=="no", ])

filteredData[filteredData$allstudents>200, ]

11. Descriptive statistics

summary(data)
##        X           minority              age           gender         
##  Min.   :  1.0   Length:463         Min.   :29.00   Length:463        
##  1st Qu.:116.5   Class :character   1st Qu.:42.00   Class :character  
##  Median :232.0   Mode  :character   Median :48.00   Mode  :character  
##  Mean   :232.0                      Mean   :48.37                     
##  3rd Qu.:347.5                      3rd Qu.:57.00                     
##  Max.   :463.0                      Max.   :73.00                     
##    credits              beauty                eval         division        
##  Length:463         Min.   :-1.4504940   Min.   :2.100   Length:463        
##  Class :character   1st Qu.:-0.6562689   1st Qu.:3.600   Class :character  
##  Mode  :character   Median :-0.0680143   Median :4.000   Mode  :character  
##                     Mean   : 0.0000001   Mean   :3.998                     
##                     3rd Qu.: 0.5456024   3rd Qu.:4.400                     
##                     Max.   : 1.9700230   Max.   :5.000                     
##     native             tenure             students       allstudents    
##  Length:463         Length:463         Min.   :  5.00   Min.   :  8.00  
##  Class :character   Class :character   1st Qu.: 15.00   1st Qu.: 19.00  
##  Mode  :character   Mode  :character   Median : 23.00   Median : 29.00  
##                                        Mean   : 36.62   Mean   : 55.18  
##                                        3rd Qu.: 40.00   3rd Qu.: 60.00  
##                                        Max.   :380.00   Max.   :581.00  
##       prof      
##  Min.   : 1.00  
##  1st Qu.:20.00  
##  Median :44.00  
##  Mean   :45.43  
##  3rd Qu.:70.50  
##  Max.   :94.00
library(moments)
library(dplyr)
## 
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
df = data.frame(age=data$age, beauty=data$beauty, eval=data$eval)
summarize(df, across(age:beauty, list(kurtosis = kurtosis, skewness = skewness)))

2.1 Выбор категоризующей переменной

В качестве категоризующего признака возьмем minority

2.2 Boxplot

data['allstudents_log'] = log(data$allstudents)
data['students_log'] = log(data$students)
dfcomp <- data %>% filter( minority == "no" | minority == "yes")

library(lattice)
bwplot(students_log ~ minority, data = dfcomp, col = c("forestgreen", "gold"), main = "students_log", xlab = "minority")

bwplot(allstudents_log ~ minority, data = dfcomp, col = c("forestgreen", "gold"), main = "allstudents_log", xlab = "minority")

bwplot(eval ~ minority, data = dfcomp, col = c("forestgreen", "gold"), main = "eval", xlab = "minority")

bwplot(beauty ~ minority, data = dfcomp, col = c("forestgreen", "gold"), main = "beauty", xlab = "minority")

bwplot(age ~ minority, data = dfcomp, col = c("forestgreen", "gold"), main = "age", xlab = "minority")

На всех выше приведённых boxplot видим, что распределения не симметричные, их медианы не совпадают, а дисперсии сильно отлличаются друг от друга.

2.3 Нормальность признаков

library(tidyr)
qqmath(~value | name, data = pivot_longer(dfcomp, c(students_log,allstudents_log)), subset = minority == "yes")

qqmath(~value | name, data = pivot_longer(dfcomp, c(students_log,allstudents_log)), groups = minority, auto.key = TRUE,
       prepanel = prepanel.qqmathline,
       panel = function(x, ...) {
          panel.qqmathline(x, ...)
          panel.qqmath(x, ...)
       })

shapiro.test(subset(dfcomp, minority == "yes")$allstudents_log)  #Критерий шапиро-уилка позволяет проверить гипотезу о том, что случайная величина имеет нормальное распределение.
## 
##  Shapiro-Wilk normality test
## 
## data:  subset(dfcomp, minority == "yes")$allstudents_log
## W = 0.89149, p-value = 3.882e-05
qqmath(~value | name, data = pivot_longer(dfcomp, c(students_log,allstudents_log)), subset = minority == "no")

shapiro.test(subset(dfcomp, minority == "no")$allstudents_log)  #Критерий шапиро-уилка позволяет проверить гипотезу о том, что случайная величина имеет нормальное распределение.
## 
##  Shapiro-Wilk normality test
## 
## data:  subset(dfcomp, minority == "no")$allstudents_log
## W = 0.93998, p-value = 1.297e-11

Критерий показал, что предположение о нормальности выбранных нами признаков отвергнуто.

2.4 t-test, критерий Манна-Уитни

t.test(allstudents_log ~ minority, data = dfcomp) #Т-тест критерий позволяет проверить гипотезу о равенстве средних по одному признаку.
## 
##  Welch Two Sample t-test
## 
## data:  allstudents_log by minority
## t = 2.8329, df = 98.711, p-value = 0.005593
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
##  0.0808376 0.4588910
## sample estimates:
##  mean in group no mean in group yes 
##          3.603930          3.334066
wilcox.test(allstudents_log ~ minority, data = dfcomp) #непараметрический аналог, работает аналогичным образом.
## 
##  Wilcoxon rank sum test with continuity correction
## 
## data:  allstudents_log by minority
## W = 14948, p-value = 0.02824
## alternative hypothesis: true location shift is not equal to 0

Здесь видим, что \(p-value<0.05\), вследствии чего гипотеза о равенстве средних отвергается. Критерий Вилкоксона менее точный, зато более устойчивый.

2.5 Критерий Колмогорова-Смирнова

dfcomp2 = subset(dfcomp, select = -c(X, gender, credits, division, native, tenure, prof) )
ks.test(dfcomp2[dfcomp$minority == "yes",6], dfcomp2[dfcomp$minority == "no", 6]) #Критерий Колмогорова-Смирнова для проверки гипотезы о равенстве функций распределения двух случайных величин.
## Warning in ks.test(dfcomp2[dfcomp$minority == "yes", 6], dfcomp2[dfcomp$minority
## == : p-value will be approximate in the presence of ties
## 
##  Two-sample Kolmogorov-Smirnov test
## 
## data:  dfcomp2[dfcomp$minority == "yes", 6] and dfcomp2[dfcomp$minority == "no", 6]
## D = 0.28074, p-value = 0.0003352
## alternative hypothesis: two-sided

Здесь видим, что \(p-value<0.05\), вследствии чего гипотеза о равенстве функций распределений отвергается.